*Stata version 15.1
*Following packages needed
net install st0085_2, from(http://www.stata-journal.com/software/sj14-2)
net install gr0002_3, from(http://www.stata-journal.com/software/sj4-3)
net install cem, from(http://fmwww.bc.edu/RePEc/bocode/c)
net install st0500, from(http://www.stata-journal.com/software/sj17-4)

cd "[directory]" // Insert directory here

********************Generate variable for analysis******************
******Generate dependent variables***
use data_municipalities.dta, clear
foreach num of numlist 2005 2009 2013 2017 {
gen diff_`num'=turnout`num'-turnout1997 // For main analysis - difference in turnout from 1997 to year
}

foreach num of numlist 2005 2009 2013 2017 {
gen diff93_`num'=turnout`num'-turnout1993 // Difference in turnout from 1993 to year
}

foreach num of numlist 2005 2009 2013 2017 {
gen diff01_`num'=turnout`num'-turnout2001 // Difference in turnout from 2001 to year
}

foreach num of numlist 1985 1989 1993 {
gen diff_`num'=turnout1997-turnout`num' // Difference in turnout from year to 1997 (for placebo analyses)
}

**Generate log-variables for level controls*
foreach num of numlist 1985 1989 1993 1997 2001 2005 2009 2013 2017 {
gen log_noeligible`num'=log(noeligible`num')
gen log_area`num'=log(area`num')
}

**Difference controls compared to 1997*
foreach num of numlist 2005 2009 2013 2017 {
gen diff_pct_youngolder`num'=pct_youngolder`num'-pct_youngolder1997 // Changes in percentage points for percent variables
gen diff_pct_immigrants`num'=pct_immigrants`num'-pct_immigrants1997
gen diff_pct_schooleduc`num'=pct_schooleduc`num'-pct_schooleduc1997
gen diff_pct_noeligible`num'=((noeligible`num'-noeligible1997)*100)/((noeligible1997)) // Changes in percent of start year
}

**Difference controls compared to 1993*
foreach num of numlist 2005 2009 2013 2017 {
gen diff93_pct_youngolder`num'=pct_youngolder`num'-pct_youngolder1993 // For robustness of alternative start year
gen diff93_pct_immigrants`num'=pct_immigrants`num'-pct_immigrants1993
gen diff93_pct_schooleduc`num'=pct_schooleduc`num'-pct_schooleduc1993
gen diff93_pct_noeligible`num'=((noeligible`num'-noeligible1993)*100)/((noeligible1993))
}

**Difference controls compared to 2001*
foreach num of numlist 2005 2009 2013 2017 {
gen diff01_pct_youngolder`num'=pct_youngolder`num'-pct_youngolder2001 // For robustness of alternative start year
gen diff01_pct_immigrants`num'=pct_immigrants`num'-pct_immigrants2001
gen diff01_pct_schooleduc`num'=pct_schooleduc`num'-pct_schooleduc2001
gen diff01_pct_noeligible`num'=((noeligible`num'-noeligible2001)*100)/((noeligible2001))
}

**Difference controls with 1997 as end year*
foreach num of numlist 1985 1989 1993 {
gen diff_pct_youngolder`num'=pct_youngolder1997-pct_youngolder`num' // For analyses of 1997 as end year (placebo analyses)
gen diff_pct_immigrants`num'=pct_immigrants1997-pct_immigrants`num'
gen diff_pct_schooleduc`num'=pct_schooleduc1997-pct_schooleduc`num'
gen diff_pct_noeligible`num'=((noeligible1997-noeligible`num')*100)/((noeligible`num'))
}

save data_analysis, replace

*********************************************************************************************
********Figure 1 - descriptive graphs for overall relationship btw size and turnout**********
*********************************************************************************************
use electiondata_oldmunicipalities, clear
foreach num of numlist 1993 1997 {
gen log_noeligible`num'=log(noeligible`num')
reg turnout`num' log_noeligible`num'
twoway (scatter turnout`num' log_noeligible`num') (lfit turnout`num' log_noeligible`num'), saving(figure1_`num', replace) ylabel(55(5)90) xlabel(7(1)13) scheme(lean1) xtitle("") ytitle("") title("`num' election") legend(off) // Elections before the reform for the graph
}

use data_municipalities.dta, clear
foreach num of numlist 2005 2009 2013 2017 {
gen log_noeligible`num'=log(noeligible`num')
reg turnout`num' log_noeligible`num'
twoway (scatter turnout`num' log_noeligible`num') (lfit turnout`num' log_noeligible`num'), saving(figure1_`num', replace) ylabel(55(5)90) xlabel(7(1)13) scheme(lean1) xtitle("") ytitle("") title("`num' election") legend(off) // Elections after the reform for the graph 
}
graph combine "figure1_1993" "figure1_1997" "figure1_2005" "figure1_2009" "figure1_2013" "figure1_2017", scheme(lean1) l1(Turnout in percent) b1(No. of eligible individuals (Ln))

graph export Figure1.tif, replace  height(6750) width(9000)

*********************************************************************************************
**************Figure 2 - descriptive graphs for treatment and control************************
*********************************************************************************************
keep munid munamalgamation turnout* 
reshape long turnout, i(munid) j(year)
bysort year: tabstat turnout, by(munamalgamation)
graph bar (mean) turnout, over(munamalgamation)  asy over(year)  ytitle("Turnout in percent")  legend(label (1 "Non-amalgamated") label(2 "Amalgamated") pos(6) row(1)) exclude0 yscale(range(60 90)) ylabels(60(5)90) scheme(Lean1)
graph export Figure2.tif, replace height(6750) width(9000)

*********************************************************************************************
*******************************Table 1 - main analyses***************************************
*********************************************************************************************
use data_analysis.dta, clear
sum munamalgamation diff_pct_immigrants2005 diff_pct_youngolder2005 diff_pct_schooleduc2005 diff_pct_noeligible2005 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

eststo q1: reg diff_2005 munamalgamation
eststo q2: reg diff_2005 munamalgamation  diff_pct_immigrants2005 diff_pct_youngolder2005 diff_pct_schooleduc2005 diff_pct_noeligible2005 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

eststo q3: reg diff_2009 munamalgamation
eststo q4: reg diff_2009 munamalgamation  diff_pct_immigrants2009 diff_pct_youngolder2009 diff_pct_schooleduc2009 diff_pct_noeligible2009 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

eststo q5: reg diff_2013 munamalgamation
eststo q6: reg diff_2013 munamalgamation  diff_pct_immigrants2013 diff_pct_youngolder2013 diff_pct_schooleduc2013 diff_pct_noeligible2013 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

eststo q7: reg diff_2017 munamalgamation, r
eststo q8: reg diff_2017 munamalgamation  diff_pct_immigrants2017 diff_pct_youngolder2017 diff_pct_schooleduc2017 diff_pct_noeligible2017 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

esttab q1 q2 q3 q4 q5 q6 q7 q8 using table1.rtf, replace varlabels(_cons "Constant" munamalgamation "Amalgamation" diff_pct_immigrants2005 "Diff immigrants pp." diff_pct_youngolder2005 "Diff young and older pp." diff_pct_schooleduc2005 ///
"Diff school education pp." diff_pct_noeligible2005 "Diff eligibles in percent" pct_immigrants1997 "Immigrants percent (1997)" pct_schooleduc1997 "School education percent (1997)" pct_youngolder1997 "Young and older percent (1997)" ///
log_noeligible1997 "Ln(eligible individuals) (1997)" log_area1997 "Ln(area (1997))" smallisland "Small island municipality") scalars("N N" "F F-test" "r2 R^2" "rmse RMSE") cells("b(star fmt(a2))" se(par)) starlevels(* 0.05 ** 0.01 *** 0.001) varwidth(50)

*********************************************************************************************
*******************************Table A1 - placebo********************************************
*********************************************************************************************
eststo qA9:  reg diff_1985 munamalgamation
eststo qA10: reg diff_1985 munamalgamation  diff_pct_immigrants1985 diff_pct_youngolder1985 diff_pct_schooleduc1985 diff_pct_noeligible1985 pct_immigrants1985 pct_youngolder1985 pct_schooleduc1985 log_noeligible1985 log_area1985 smallisland

eststo qA11: reg diff_1989 munamalgamation, r
eststo qA12: reg diff_1989 munamalgamation  diff_pct_immigrants1989 diff_pct_youngolder1989 diff_pct_schooleduc1989 diff_pct_noeligible1989 pct_immigrants1989 pct_youngolder1989 pct_schooleduc1989 log_noeligible1989 log_area1989 smallisland 

eststo qA13: reg diff_1993 munamalgamation
eststo qA14: reg diff_1993 munamalgamation  diff_pct_immigrants1993 diff_pct_youngolder1993 diff_pct_schooleduc1993 diff_pct_noeligible1993 pct_immigrants1993 pct_youngolder1993 pct_schooleduc1993 log_noeligible1993 log_area1993 smallisland

esttab qA9 qA10 qA11 qA12 qA13 qA14 using tableA1.rtf, replace scalars("N N" "F F-test" "r2 R^2" "rmse RMSE") cells("b(star fmt(a2))" se(par)) starlevels(* 0.05 ** 0.01 *** 0.001) varwidth(50)

*********************************************************************************************
*******************************Table A2 - matching********************************************
*********************************************************************************************
cem noeligible1997 (5000 10000 20000 30000 40000 70000 100000 200000), treatment(munamalgamation) 
eststo qA15: reg diff_2005 munamalgamation diff_pct_immigrants2005 diff_pct_youngolder2005 diff_pct_schooleduc2005 diff_pct_noeligible2005 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland [aw=cem_weights]
eststo qA16: reg diff_2009 munamalgamation diff_pct_immigrants2009 diff_pct_youngolder2009 diff_pct_schooleduc2009 diff_pct_noeligible2009 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland [aw=cem_weights]
eststo qA17: reg diff_2013 munamalgamation diff_pct_immigrants2013 diff_pct_youngolder2013 diff_pct_schooleduc2013 diff_pct_noeligible2013 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland [aw=cem_weights]
eststo qA18: reg diff_2017 munamalgamation diff_pct_immigrants2017 diff_pct_youngolder2017 diff_pct_schooleduc2017 diff_pct_noeligible2017 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland [aw=cem_weights]

esttab qA15 qA16 qA17 qA18 using tableA2.rtf, replace scalars("N N" "F F-test" "r2 R^2" "rmse RMSE") cells("b(star fmt(a2))" se(par)) starlevels(* 0.05 ** 0.01 *** 0.001) varwidth(50)

*********************************************************************************************
****************************Table A3 - diff-in-diff******************************************
*********************************************************************************************
preserve
keep munid munamalgamation turnout19* turnout20* pct* log_area* log_noelig* smallisland 
order munid turnout* munamalgamation
reshape long turnout pct_immigrants pct_youngolder pct_schooleduc log_noeligible log_area, i(munid) j(year) // Data in long format for these analyses

eststo qA19: reg turnout b1997.year##i.munamalgamation , cluster(munid)
eststo qA20: reg turnout b1997.year##i.munamalgamation  pct_immigrants pct_youngolder pct_schooleduc log_noeligible log_area smallisland, cluster(munid)
eststo qA21: reg turnout b1997.year##i.munamalgamation  b1997.year##c.pct_immigrants b1997.year##c.pct_youngolder b1997.year##c.pct_schooleduc b1997.year##c.log_noeligible b1997.year##c.log_area b1997.year##c.smallisland, cluster(munid)

esttab qA19 qA20 qA21 using tableA3.rtf, replace varlabels(_cons "Constant" pct_immigrants "Immigrants percent" pct_schooleduc "School education percent" ///
pct_youngolder "Young and older percent" log_noeligible "Ln(eligible individuals)"  log_area "Ln(Area)" smallisland "Small island municipality") ///
scalars("N N" "F F-test" "r2 R^2" "rmse RMSE") cells("b(star fmt(a2))" se(par)) starlevels(* 0.05 ** 0.01 *** 0.001) varwidth(50)
restore

*********************************************************************************************
*******************************Table A4 - alternative start year*****************************
*********************************************************************************************
eststo qA22: reg diff93_2005 munamalgamation, r
eststo qA23: reg diff93_2005 munamalgamation  diff93_pct_immigrants2005 diff93_pct_youngolder2005 diff93_pct_schooleduc2005 diff93_pct_noeligible2005 pct_immigrants1993 pct_youngolder1993 pct_schooleduc1993 log_noeligible1993 log_area1993 smallisland

eststo qA24: reg diff93_2009 munamalgamation
eststo qA25: reg diff93_2009 munamalgamation  diff93_pct_immigrants2009 diff93_pct_youngolder2009 diff93_pct_schooleduc2009 diff93_pct_noeligible2009 pct_immigrants1993 pct_youngolder1993 pct_schooleduc1993 log_noeligible1993 log_area1993 smallisland

eststo qA26: reg diff93_2013 munamalgamation, r
eststo qA27: reg diff93_2013 munamalgamation  diff93_pct_immigrants2013 diff93_pct_youngolder2013 diff93_pct_schooleduc2013 diff93_pct_noeligible2013 pct_immigrants1993 pct_youngolder1993 pct_schooleduc1993 log_noeligible1993 log_area1993 smallisland

eststo qA28: reg diff93_2017 munamalgamation, r
eststo qA29: reg diff93_2017 munamalgamation  diff93_pct_immigrants2017 diff93_pct_youngolder2017 diff93_pct_schooleduc2017 diff93_pct_noeligible2017 pct_immigrants1993 pct_youngolder1993 pct_schooleduc1993 log_noeligible1993 log_area1993 smallisland

esttab qA22 qA23 qA24 qA25 qA26 qA27 qA28 qA29 using tableA4.rtf, replace scalars("N N" "F F-test" "r2 R^2" "rmse RMSE") cells("b(star fmt(a2))" se(par)) starlevels(* 0.05 ** 0.01 *** 0.001) varwidth(50)

*********************************************************************************************
*************************Table A6 - amalgamation degree**************************************
*********************************************************************************************
eststo qA30: reg diff_2005 herfindahl
eststo qA31: reg diff_2005 herfindahl  diff_pct_immigrants2005 diff_pct_youngolder2005 diff_pct_schooleduc2005 diff_pct_noeligible2005 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

eststo qA32: reg diff_2009 herfindahl
eststo qA33: reg diff_2009 herfindahl  diff_pct_immigrants2009 diff_pct_youngolder2009 diff_pct_schooleduc2009 diff_pct_noeligible2009 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

eststo qA34: reg diff_2013 herfindahl
eststo qA35: reg diff_2013 herfindahl  diff_pct_immigrants2013 diff_pct_youngolder2013 diff_pct_schooleduc2013 diff_pct_noeligible2013 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

eststo qA36: reg diff_2017 herfindahl
eststo qA37: reg diff_2017 herfindahl  diff_pct_immigrants2017 diff_pct_youngolder2017 diff_pct_schooleduc2017 diff_pct_noeligible2017 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

esttab qA30 qA31 qA32 qA33 qA34 qA35 qA36 qA37 using tableA6.rtf, replace scalars("N N" "F F-test" "r2 R^2" "rmse RMSE") cells("b(star fmt(a2))" se(par)) starlevels(* 0.05 ** 0.01 *** 0.001) varwidth(50)

*********************************************************************************************
*************************Table A7 - amalgamation degree squared******************************
*********************************************************************************************
eststo qA38: reg diff_2005 c.herfindahl##c.herfindahl
eststo qA39: reg diff_2005 c.herfindahl##c.herfindahl  diff_pct_immigrants2005 diff_pct_youngolder2005 diff_pct_schooleduc2005 diff_pct_noeligible2005 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

eststo qA40: reg diff_2009 c.herfindahl##c.herfindahl
eststo qA41: reg diff_2009 c.herfindahl##c.herfindahl  diff_pct_immigrants2009 diff_pct_youngolder2009 diff_pct_schooleduc2009 diff_pct_noeligible2009 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

eststo qA42: reg diff_2013 c.herfindahl##c.herfindahl
eststo qA43: reg diff_2013 c.herfindahl##c.herfindahl  diff_pct_immigrants2013 diff_pct_youngolder2013 diff_pct_schooleduc2013 diff_pct_noeligible2013 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

eststo qA44: reg diff_2017 c.herfindahl##c.herfindahl
eststo qA45: reg diff_2017 c.herfindahl##c.herfindahl  diff_pct_immigrants2017 diff_pct_youngolder2017 diff_pct_schooleduc2017 diff_pct_noeligible2017 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland

esttab qA38 qA39 qA40 qA41 qA42 qA43 qA44 qA45 using tableA7.rtf, replace scalars("N N" "F F-test" "r2 R^2" "rmse RMSE") cells("b(star fmt(a2))" se(par)) starlevels(* 0.05 ** 0.01 *** 0.001) varwidth(50)

*********************************************************************************************
*************************Table A8 - amalgamation degree alternative measure******************
*********************************************************************************************
recode nooldmunicipalities (1/1=1) (2=2) (3/4=3) (5/7=4), gen(nooldmunicipalities_dummies)
eststo qA46: reg diff_2005 i.nooldmunicipalities_d
eststo qA47: reg diff_2005 i.nooldmunicipalities_d  diff_pct_immigrants2005 diff_pct_youngolder2005 diff_pct_schooleduc2005 diff_pct_noeligible2005 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland
testparm i.nooldmunicipalities_dummies

eststo qA48: reg diff_2009 i.nooldmunicipalities_d
eststo qA49: reg diff_2009 i.nooldmunicipalities_d  diff_pct_immigrants2009 diff_pct_youngolder2009 diff_pct_schooleduc2009 diff_pct_noeligible2009 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland
testparm i.nooldmunicipalities_dummies

eststo qA50: reg diff_2013 i.nooldmunicipalities_d
eststo qA51: reg diff_2013 i.nooldmunicipalities_d  diff_pct_immigrants2013 diff_pct_youngolder2013 diff_pct_schooleduc2013 diff_pct_noeligible2013 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland
testparm i.nooldmunicipalities_dummies

eststo qA52: reg diff_2017 i.nooldmunicipalities_d, r
eststo qA53: reg diff_2017 i.nooldmunicipalities_d diff_pct_immigrants2017 diff_pct_youngolder2017 diff_pct_schooleduc2017 diff_pct_noeligible2017 pct_immigrants1997 pct_youngolder1997 pct_schooleduc1997 log_noeligible1997 log_area1997 smallisland
testparm i.nooldmunicipalities_dummies

esttab qA46 qA47 qA48 qA49 qA50 qA51 qA52 qA53 using tableA8.rtf, replace scalars("N N" "F F-test" "r2 R^2" "rmse RMSE") cells("b(star fmt(a2))" se(par)) starlevels(* 0.05 ** 0.01 *** 0.001) varwidth(50)

*********************************************************************************************
***************************Table A5 - synthetic control**************************************
*********************************************************************************************
***Analysis in long format**
keep munid munamalgamation turnout19* turnout20* pct* log_area* log_noelig* smallisland
order munid turnout* munamalgamation
reshape long turnout pct_immigrants pct_youngolder pct_schooleduc log_noeligible log_area, i(munid) j(year)
drop if year==2001

tsset munid year, delta(4)

**Using synth_runner*
gen treatment=munamalgamation
replace treatment=0 if year<2005
synth_runner turnout pct_immigrants pct_youngolder pct_schooleduc log_noeligible log_area smallisland turnout(1985) turnout(1997), d(treatment)

*Main robustness analysis*
foreach num of numlist 190 210 219 230 240 250 259 260 265 270 306 316 320 326 330 336 340 350 360 370 376 390 400 410 420 430 440 ///
450 479 480 482 492 510 530 540 550 561 573 575 580 615 621 630 657 661 665 671 706 707 710 730 740 746 756 760 766 779 787 791 810 813 820 840 846 849 851 860 {
synth turnout pct_immigrants pct_youngolder pct_schooleduc log_noeligible log_area smallisland turnout(1985) turnout(1997), trunit(`num') ///
counit(101 147 151 153 155 157 159 161 163 165 167 169 173 175 183 185 187 201 217 223 253 269 329 461 563 607 727 741 751 773 825) trperiod(2005) keep(res_`num', replace) 

preserve
use res_`num', clear
mat rmspe=e(RMSPE)
svmat rmspe, names(rmspe)
replace rmspe1=rmspe1[_n-1] if missing(rmspe1)
drop _Co_Number _W_Weight
gen munid=`num'
save res_`num', replace
restore
}

use res_190, clear
foreach num of numlist 210 219 230 240 250 259 260 265 270 306 316 320 326 330 336 340 350 360 370 376 390 400 410 420 430 440 ///
450 479 480 482 492 510 530 540 550 561 573 575 580 615 621 630 657 661 665 671 706 707 710 730 740 746 756 760 766 779 787 791 810 813 820 840 846 849 851 860 {
append using res_`num'
}

gen diff_Y=_Y_treated -_Y_synthetic 

bysort _time: sum _Y_synthetic _Y_treated diff_Y // Table A5
